Background

Single cells were clustered into 9 clusters based on the expression of proteins that were not differentially expressed between drug treatments. These clusters appear to have some relation to cell-state. And interesting question is if some phospho-protein responses are cluster dependent, as this would suggest a cell-state dependent drug response that is potentially relevant to drug resistance.

In this notebook we systematically explore this. To this end, we regress the expression of each antibody in each cell to the drug-treatment and cellstate cluster, and include an interaction term, i.e. we fit the following model for each antibody.

\[expression_{i,x} ~ \beta_0 + \beta_{t}*treatment + \beta_{c}*cellstate + \beta_{i}*treatment * cellstate\]

Antibodies for which there are significant interaction term show signs of cell-state dependent drug response.

We will use the z-score normalized data as input because this allows for visualization of all antibodies on the same scale. Note that this does not the p-values for significance of the interaction term (which is shown at the end of the notebook).

Our input data looks as follows:

dat_tmm <- read_tsv(here("data", "processed", "scIDseq-vanEijl-tmm.tsv")) %>% 
  select(sample_id, treatment, ab_name, ab_type, ab_count_tmm) %>% 
  filter(treatment != "No_EGF") %>% 
  mutate(treatment = factor(treatment, levels = c("EGF", "ip70S6K_EGF", "iRSK_EGF")))

dat_zscore <- read_csv(
  here("data", "processed" ,"TMM_normalised_z_transformed_concencus_clusters.csv")
  ) %>% 
  select(-X1, -plate_number, -treatment) %>% 
  mutate(cluster = as_factor(cluster)) %>% 
  pivot_longer(!c(cluster, sample_id), names_to = "ab_name", values_to = "zscore")


dat <- left_join(dat_tmm, dat_zscore)

rm(dat_tmm, dat_zscore)
head(dat)
length(unique(dat$sample_id))
[1] 546
plot_ab <- function(ab) {
  ggplot(
    filter(dat, ab_name == ab),
    aes(x = cluster, y = zscore, fill = treatment)
  ) +
    geom_boxplot() +
    #expand_limits(y = 0) +
    my_theme +
    scale_fill_manual(values = col_lst) +
    labs(title = ab, x = "Cell state cluster", y = "z-score normalized ab-expression")
}

plot_ab("RIBOSOMAL_S6_P")

Analysis

Data exploration

Let’s first briefly explore the clusters. How many cells are there in each cluster, and which treatment did they receive?

dat %>% 
  ungroup() %>% 
  select(sample_id, treatment, cluster) %>% 
  group_by(treatment, cluster) %>% 
  distinct() %>% 
  summarise(N = n()) %>% 
  pivot_wider(names_from = treatment, values_from = N) %>%
  mutate(total = EGF + ip70S6K_EGF + iRSK_EGF)

Or visually.

dat %>% 
  select(sample_id, treatment, cluster) %>% 
  distinct() %>% 
  ggplot(aes(x = cluster, fill = treatment)) +
  geom_bar() +
  scale_fill_manual(values = col_lst) +
  my_theme

We expect the treatments to be balanced over the clusters, as the duration of the treatment is too short to affect the cell state. This seems to be reasonably the case, but a formal enrichment analysis might be useful here.

Model fitting

Now let’s fit our model with interaction term. In addition, we obtain the “fit” for each treatment-cluster combination and the corresponding the 95% confidence interval. For comparison, we’ll also fit a model without an interaction term.

As reference levels, we select EGF only treatment and cluster 1. Later we’ll explore the effect of using different clusters as reference.

get_fit <- function(ab, df, formula){
  lm(formula, filter(df, ab_name == ab))
}


formula_simple <- zscore ~  0 + cluster + treatment
formula_interaction <- zscore ~ 0 + cluster * treatment

ab_name_lst <- as.list(unique(dat$ab_name))
names(ab_name_lst) <- ab_name_lst

fit_lst_simple <- purrr::map(ab_name_lst, get_fit, dat, formula_simple)
fit_lst_interaction <- purrr::map(ab_name_lst, get_fit, dat, formula_interaction)

model_simple <- 
  purrr::map(fit_lst_simple, tidy) %>% 
  bind_rows(.id = "ab_name") %>% 
  mutate(fdr = p.adjust(p.value, method = "BH")) %>% 
  mutate(term_type = case_when(
    str_detect(term, ":") ~ "interaction",
    str_detect(term, "(Intercept)") ~ "intercept",
    str_detect(term, "treatment") ~ "treatment",
    TRUE ~ "cluster"
  )) 

model_interaction <- 
  purrr::map(fit_lst_interaction, tidy) %>% 
  bind_rows(.id = "ab_name") %>% 
  mutate(fdr = p.adjust(p.value, method = "BH")) %>% 
  mutate(term_type = case_when(
    str_detect(term, ":") ~ "interaction",
    str_detect(term, "(Intercept)") ~ "intercept",
    str_detect(term, "treatment") ~ "treatment",
    TRUE ~ "cluster"
  )) 

write_tsv(model_simple, here("results", "model-exp-lm"))
write_tsv(model_interaction, here("results", "model-exp-interaction.tsv"))

Two-way Anova

We’ll first perform a two-way ANOVA as this is the simplest to interpret. This will give us significance levels per antibody, but not for individual interactions terms. The following antibodies have a significant interaction term (which means that the zscores can be significantly better described by a model with an interaction term compared to one without.)

model_anova <- purrr::map(
  ab_name_lst,
  function(x){tidy(aov(zscore ~ treatment * cluster, data = filter(dat, ab_name == x)))}
  ) %>% 
  bind_rows(.id = "ab_name")  

signif_interactions_anova <- 
  model_anova %>% 
  filter(term == "treatment:cluster") %>% 
  mutate(fdr = p.adjust(p.value, method = "BH")) %>% 
  arrange(fdr) %>% 
  filter(fdr < 0.05) 
  
select(signif_interactions_anova, ab_name, statistic, p.value, fdr)
write_tsv(signif_interactions_anova, here("results",  "treatment_cellstate_interactions_anova.tsv"))

Regression with interaction term

To get p-values for each individual interaction term we .

Note that the final model is mathematically equivalent to the two-way ANOVA in that model predictions/fits are identical. Let’s consider the significant interactions. A table with all significant interactions is written to results/scIDseq/treatment_cellstate_interactions.tsv, and a table with all coefficients of the linear and interaction model to results/scIDseq/model-exp-lm.tsv and results/scIDseq/model-exp-interaction.tsv, respectively.

signif_interactions <-
  model_interaction %>% 
  filter(term_type == "interaction" & p.value < 0.05) %>%
  arrange(fdr) %>% 
  separate(term, into = c("cluster", "treatment"), sep = ":", remove = FALSE) %>% 
  mutate(
    treatment = stringr::str_replace(treatment, "treatment", ""),
    cluster = stringr::str_replace(cluster, "cluster", "")
  )  

write_tsv(signif_interactions, here("results", "treatment_cellstate_interactions.tsv"))
filter(signif_interactions, cluster == "5")

There are significant interactions. What are the stongest ones?

head(signif_interactions, n = 10)

Cluster 9 seems to show the most. Indeed, most are between p70S6K inhibition and cluster 9.

group_by(signif_interactions, cluster, treatment) %>% 
  summarise(N = n()) %>% 
  arrange(desc(N))

However, it should be noted that this is based on just 13 cells.

Visualizations

prediction_input <- distinct(select(dat, treatment, cluster)) %>% 
  arrange(cluster, treatment)

get_predictions <- function(fit){
  bind_cols(
    prediction_input,
    as_tibble(predict(fit, prediction_input, interval = "confidence")))
}

predictions_interaction <- 
  fit_lst_interaction %>% 
  purrr::map(get_predictions) %>% 
  bind_rows(.id = "ab_name")

predictions_simple <- 
  fit_lst_simple %>% 
  purrr::map(get_predictions) %>% 
  bind_rows(.id = "ab_name")


predictions <- bind_rows(list(
  "interaction" = predictions_interaction,
  "simple" = predictions_simple),
  .id = "model_type"
  ) 
rm(predictions_interaction, predictions_simple, prediction_input)

Let’s visualize the expression patterns of some potentially interesting antibodies. In the plots below the model predictions for each treatment-cluster pair are indicated. For comparison, the gray points indicate the estimates using a model without an interaction term.

plot_prediction <- function(ab){
  ggplot(filter(predictions, ab_name == ab & model_type == "interaction"), 
    aes(y = fit, x = cluster, color = treatment, group=treatment)) +
    geom_point(
      data = filter(predictions, ab_name == ab & model_type == "simple"),
      color = "lightgray", size = 4, position = position_dodge(width = 0.5)) +
    geom_errorbar(
      data = filter(predictions, ab_name == ab & model_type == "simple"),
      aes(y = fit, ymin = lwr, ymax = upr, group = treatment),
      width = .5, position = position_dodge(width = 0.5), color = "lightgray"
    ) +
    geom_point(size = 4, position = position_dodge(width = 0.5)) +
    geom_errorbar(
      aes(y = fit, ymin = lwr, ymax = upr, color = treatment),
      width = .5, position = position_dodge(width = 0.5)
    ) +
  scale_color_manual(values = col_lst) + 
  labs(title = ab, x = "Cell state cluster", y = "z-score (Model estimate)") +
  my_theme
}

RSKi treatment and Cluster 5

Cluster 5 has high expression of both H3_P and CYCLIN_B1 (and to a lesser extend, H2A_P), which is affected by RSKi. For H3_P and H2A_P the upregulation in cluster 5 is repressed, whereas for CYCLIN_B1 it is extra pronounced

plot_prediction("H3_P")

plot_ab("H3_P")

plot_prediction("H2A_P")

plot_prediction("CYCLIN_B1")

p70i treatment and cluster 9

A large number of antibodies show a similar trend, where they are upregulated upon p70S6K inhibition in cluster 9. These include CYCLIN_E, FAK, SCR_P, JAK1_P, STAT5_P, and phospho-NFkB/p65 all show a SRC_P also shows this trent in cluster 6. Surprisingly, AKT2 also shows this profile, since the response of a total protein is not expected on the timescale of the inhibitor treatment (2 hours?)

plot_prediction("CYCLIN_E")

plot_prediction("FAK")

plot_prediction("SRC_P")

plot_prediction("JAK1_P")

plot_prediction("STAT5_P")

plot_prediction("EGFR_P_Y1045")

plot_prediction("AKT2")

Other interesting antibodies

Finally, RIBOSOMAL_S6_P also has a significant interaction term, which is not so easy to pin to one cluster. In clusters 2, it is unexpectedly unresponsive to p70 inhbition, whereas in cluster 6 it is unresponsive to RSK inhibition.

plot_prediction("RIBOSOMAL_S6_P")

Overview of all term

We’ll also create an overview of all interaction terms. Signifance is indicated by a black circle around the dot.

ordered_ab_lst <- read_csv(
  here("data", "for-visualization", "ordered_antibody_list.csv"), 
   col_names = FALSE)$X1

hc_treats <- filter(model_interaction, term_type == "interaction") %>% 
  separate(term, into = c("cluster", "treatment"), sep = ":", remove = FALSE) %>% 
  mutate(log10p = -log10(p.value)) %>% 
  select(ab_name, term, log10p) %>% 
  pivot_wider(names_from = term, values_from = log10p) %>% 
  column_to_rownames("ab_name") %>% 
  #as.matrix %>% 
  dist() %>% 
  hclust()

plt <- model_interaction %>% 
  filter(term_type == "interaction") %>%
  mutate(significant = p.value < 0.05) %>%
  separate(term, into = c("cluster", "treatment"), sep = ":", remove = FALSE) %>%
  mutate(
    treatment = stringr::str_replace(treatment, "treatment", ""),
    cluster = stringr::str_replace(cluster, "cluster", "")
  )  %>% 
  mutate(
    ab_name = factor(ab_name, levels = ordered_ab_lst)
    ) %>%
  ggplot(aes(x = ab_name, y = cluster)) +
  geom_point(aes(size = -log10(p.value), fill = estimate, color = significant), shape = 21) +
  scale_size(range = c(0,5)) +
  scale_color_manual(breaks = c(TRUE, FALSE), values = c("black", "white")) +
  colorspace::scale_fill_continuous_diverging("Blue-Red 2") +
  cowplot::theme_minimal_grid(font_size = 9, font_family = "Helvetica") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
  labs(x = "", y = "Cell state") +
  guides(color = FALSE) +
  facet_wrap(~treatment, nrow = 2)
plt

ggsave(here("figures",  "interaction-terms-ref=C1.pdf"), plt, width = 10, height = 4)

Using different reference clusters

The plot above does not show any interaction terms for cluster 1 because we used cluster 1 as the reference in the plot above. However, this is a somewhat arbitrary choice. Inspection of some of the individual antibodies showed interesting interactions in Clusters, 2, 5, 6, and 9, so we don’t want to use these as reference.

One way to select a reference cluster is to use the most “average” one. A rough way to determine this is looking at the mean of the absolute z-scores of all antibodies.

dat %>% 
  group_by(cluster) %>% 
  summarise(
    mean_zscore_abs = mean(abs(zscore)),
    mean_zscore_sq = mean(zscore*zscore)
    ) %>% 
  arrange(mean_zscore_sq)

This would suggest using cluster 3 as a reference.

dat_ref3 <- mutate(dat, cluster = forcats::fct_relevel(cluster, "3"))
fit_lst_zscore_ref3 <- purrr::map(ab_name_lst, get_fit, dat_ref3, formula_interaction)

model_ref3 <-
  purrr::map(fit_lst_zscore_ref3, function(x) broom::tidy(x)) %>%
  bind_rows(.id = "ab_name") %>%
  mutate(fdr = p.adjust(p.value, method = "BH")) %>%
  mutate(term_type = case_when(
    str_detect(term, ":") ~ "interaction",
    str_detect(term, "(Intercept)") ~ "intercept",
    str_detect(term, "treatment") ~ "treatment",
    TRUE ~ "cluster"
  ))

signif_interactions_ref3 <-
  model_ref3 %>%
  filter(term_type == "interaction" & p.value < 0.05) %>%
  arrange(fdr) %>%
  separate(term, into = c("cluster", "treatment"), sep = ":", remove = FALSE) %>%
  mutate(
    treatment = stringr::str_replace(treatment, "treatment", ""),
    cluster = stringr::str_replace(cluster, "cluster", "")
  )

plt <- 
  model_ref3 %>%
  filter(term_type == "interaction") %>%
  mutate(significant = p.value < 0.05) %>%
  separate(term, into = c("cluster", "treatment"), sep = ":", remove = FALSE) %>%
  mutate(
    treatment = stringr::str_replace(treatment, "treatment", ""),
    cluster = stringr::str_replace(cluster, "cluster", "")
  )  %>%
  mutate(ab_name = factor(ab_name, levels = ordered_ab_lst)) %>% 
  # mutate(cluster = factor(cluster, levels = hc_clust$labels[hc_clust$order])) %>% 
  ggplot(aes(x = ab_name, y = cluster)) +  
  geom_point(aes(size = -log10(p.value), fill = estimate, color = significant), shape = 21) +
  scale_size(range = c(0,5)) +
  scale_color_manual(breaks = c(TRUE, FALSE), values = c("black", "white")) +
  colorspace::scale_fill_continuous_diverging("Blue-Red 2") +
  cowplot::theme_minimal_grid(font_size = 9, font_family = "Helvetica") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
  labs(x = "", y = "Cell state") +
  guides(color = FALSE) +
  facet_wrap(~treatment, nrow = 2)
plt

ggsave(here("figures",  "interaction-terms-ref=C3.pdf"), plt, width = 10, height = 4)

Of note, when using cluster 3 as reference, the interaction terms relating to RPS6 are no longer significant. (When using cluster 1 also not after multiple testing correction).

How much difference does the reference cluster make? Let’s look at the correlation between estimate effect sizes for the interaction term using Cluster 1 and cluster 3 as reference.

tmp <- filter(model_interaction, term_type == "interaction") %>% 
  full_join(
    filter(model_ref3, term_type == "interaction"), 
    by = c("ab_name", "term"), suffix = c(".ref1", ".ref3")
    ) %>% 
  separate(term, into = c("cluster", "treatement"), sep = ":")


ggplot(tmp, aes(x = estimate.ref1, y=estimate.ref3)) +
  geom_point() +
  my_theme +
  labs(
    x = "Cluster 1 as reference", y = "Cluster 3 as reference",
    title = "Effect sizes of interaction term"
    )

The plot above shows that these are highly correlated.

Using altenative inputs

We have used the TMM normalized data to perform the interaction term analysis. Here we just check that using the TMM-data gives identical p-values for the two-way anova and the regression model.

formula_tmm <- ab_count_tmm ~ cluster * treatment
fit_lst_tmm <- purrr::map(ab_name_lst, get_fit, dat, formula_tmm)

model_tmm <-
  purrr::map(fit_lst_tmm, function(x) broom::tidy(x)) %>%
  bind_rows(.id = "ab_name") %>%
  mutate(fdr = p.adjust(p.value, method = "BH")) %>%
  mutate(term_type = case_when(
    str_detect(term, ":") ~ "interaction",
    str_detect(term, "(Intercept)") ~ "intercept",
    str_detect(term, "treatment") ~ "treatment",
    TRUE ~ "cluster"
  ))

signif_interactions_tmm <-
  model_tmm %>%
  filter(term_type == "interaction" & p.value < 0.05) %>%
  arrange(fdr) %>%
  separate(term, into = c("cluster", "treatment"), sep = ":", remove = FALSE) %>%
  mutate(
    treatment = stringr::str_replace(treatment, "treatment", ""),
    cluster = stringr::str_replace(cluster, "cluster", "")
  )

sum(signif_interactions_tmm$p.value - signif_interactions$p.value > 1e-3)
[1] 0
rm()
plot_correlation <- function(df, ab1, ab2, c, count_col = 'zscore'){
  df %>%
    select(sample_id, ab_name, treatment, cluster, !!sym(count_col)) %>%
    pivot_wider(
      names_from = ab_name,
      values_from = !!sym(count_col)
    ) %>% 
    mutate(in_cluster = if_else(cluster == c, TRUE, FALSE)) %>% 
    ggplot(aes(x = !!sym(ab1), y = !!sym(ab2), col = in_cluster)) +
    geom_point() +
    facet_wrap(~treatment) +
    labs(x = ab1, y = ab2, title = str_c("Cluster ", c), subtitle = count_col) +
    geom_smooth(method = "lm") +
    scale_color_manual(values = c("gray", "red"), name = str_c("In cluster ", c))
}

plot_correlation(dat,"MAPK_APK2_P",  "CYCLIN_B1", "5") 

plot_correlation(dat,"MAPK_APK2_P",  "CYCLIN_B1", "9") 

# plot_correlation(dat,"MAPK_APK2_P",  "CYCLIN_B1", "2") 
plot_correlation(dat,"MKK3_6_P",  "CYCLIN_B1", "5") 

plot_correlation(dat, "RB_P", "H2A_P", "5") 

# plot_correlation(dat, "EGFR_P_Y1045", "ERK1_2_P")
plot_correlation(dat,"MAPK_APK2_P",  "CYCLIN_B1", "5") + facet_wrap(~treatment, ncol=1)

#plot_correlation(dat, "CJUN_P", "MAPK_APK2_P")  + facet_wrap(~treatment, scales = "free")
plot_correlation(dat, "ERK1_2_P", "MAPK_APK2_P", "5") + facet_wrap(~treatment, ncol=1)

Session info

sessioninfo::session_info()
─ Session info ───────────────────────────────────────────────────────────────
 setting  value                       
 version  R version 4.0.3 (2020-10-10)
 os       macOS Mojave 10.14.6        
 system   x86_64, darwin17.0          
 ui       X11                         
 language (EN)                        
 collate  en_US.UTF-8                 
 ctype    en_US.UTF-8                 
 tz       Europe/Amsterdam            
 date     2021-09-16                  

─ Packages ───────────────────────────────────────────────────────────────────
 package     * version    date       lib source        
 assertthat    0.2.1      2019-03-21 [1] CRAN (R 4.0.0)
 backports     1.2.1      2020-12-09 [1] CRAN (R 4.0.2)
 bookdown      0.21       2020-10-13 [1] CRAN (R 4.0.3)
 broom       * 0.7.3      2020-12-16 [1] CRAN (R 4.0.2)
 cellranger    1.1.0      2016-07-27 [1] CRAN (R 4.0.0)
 class         7.3-18     2021-01-24 [1] CRAN (R 4.0.2)
 cli           2.2.0      2020-11-20 [1] CRAN (R 4.0.2)
 codetools     0.2-18     2020-11-04 [1] CRAN (R 4.0.2)
 colorspace    2.0-0      2020-11-11 [1] CRAN (R 4.0.2)
 cowplot       1.1.1      2020-12-30 [1] CRAN (R 4.0.2)
 crayon        1.3.4      2017-09-16 [1] CRAN (R 4.0.0)
 DBI           1.1.1      2021-01-15 [1] CRAN (R 4.0.2)
 dbplyr        2.0.0      2020-11-03 [1] CRAN (R 4.0.2)
 dials       * 0.0.9      2020-09-16 [1] CRAN (R 4.0.2)
 DiceDesign    1.8-1      2019-07-31 [1] CRAN (R 4.0.0)
 digest        0.6.27     2020-10-24 [1] CRAN (R 4.0.2)
 dplyr       * 1.0.3      2021-01-15 [1] CRAN (R 4.0.2)
 ellipsis      0.3.1      2020-05-15 [1] CRAN (R 4.0.0)
 evaluate      0.14       2019-05-28 [1] CRAN (R 4.0.0)
 fansi         0.4.2      2021-01-15 [1] CRAN (R 4.0.2)
 farver        2.0.3      2020-01-16 [1] CRAN (R 4.0.0)
 forcats     * 0.5.0      2020-03-01 [1] CRAN (R 4.0.0)
 foreach       1.5.1      2020-10-15 [1] CRAN (R 4.0.2)
 fs            1.5.0      2020-07-31 [1] CRAN (R 4.0.2)
 furrr         0.2.1      2020-10-21 [1] CRAN (R 4.0.2)
 future        1.21.0     2020-12-10 [1] CRAN (R 4.0.2)
 generics      0.1.0      2020-10-31 [1] CRAN (R 4.0.2)
 ggplot2     * 3.3.3      2020-12-30 [1] CRAN (R 4.0.2)
 globals       0.14.0     2020-11-22 [1] CRAN (R 4.0.2)
 glue          1.4.2      2020-08-27 [1] CRAN (R 4.0.2)
 gower         0.2.2      2020-06-23 [1] CRAN (R 4.0.2)
 GPfit         1.0-8      2019-02-08 [1] CRAN (R 4.0.0)
 gtable        0.3.0      2019-03-25 [1] CRAN (R 4.0.0)
 haven         2.3.1      2020-06-01 [1] CRAN (R 4.0.0)
 here        * 1.0.1      2020-12-13 [1] CRAN (R 4.0.2)
 hms           1.0.0      2021-01-13 [1] CRAN (R 4.0.2)
 htmltools     0.5.1.1    2021-01-22 [1] CRAN (R 4.0.2)
 httr          1.4.2      2020-07-20 [1] CRAN (R 4.0.2)
 infer       * 0.5.4      2021-01-13 [1] CRAN (R 4.0.2)
 ipred         0.9-9      2019-04-28 [1] CRAN (R 4.0.0)
 iterators     1.0.13     2020-10-15 [1] CRAN (R 4.0.2)
 jsonlite      1.7.2      2020-12-09 [1] CRAN (R 4.0.2)
 knitr       * 1.30       2020-09-22 [1] CRAN (R 4.0.2)
 labeling      0.4.2      2020-10-20 [1] CRAN (R 4.0.2)
 lattice       0.20-41    2020-04-02 [1] CRAN (R 4.0.3)
 lava          1.6.8.1    2020-11-04 [1] CRAN (R 4.0.2)
 lhs           1.1.1      2020-10-05 [1] CRAN (R 4.0.2)
 lifecycle     0.2.0      2020-03-06 [1] CRAN (R 4.0.0)
 listenv       0.8.0      2019-12-05 [1] CRAN (R 4.0.0)
 lubridate     1.7.9.2    2020-11-13 [1] CRAN (R 4.0.2)
 magrittr      2.0.1      2020-11-17 [1] CRAN (R 4.0.2)
 MASS          7.3-53     2020-09-09 [1] CRAN (R 4.0.3)
 Matrix        1.3-2      2021-01-06 [1] CRAN (R 4.0.2)
 mgcv          1.8-33     2020-08-27 [1] CRAN (R 4.0.3)
 modeldata   * 0.1.0      2020-10-22 [1] CRAN (R 4.0.2)
 modelr        0.1.8      2020-05-19 [1] CRAN (R 4.0.0)
 munsell       0.5.0      2018-06-12 [1] CRAN (R 4.0.0)
 nlme          3.1-151    2020-12-10 [1] CRAN (R 4.0.2)
 nnet          7.3-15     2021-01-24 [1] CRAN (R 4.0.2)
 parallelly    1.23.0     2021-01-04 [1] CRAN (R 4.0.2)
 parsnip     * 0.1.5      2021-01-19 [1] CRAN (R 4.0.2)
 pillar        1.4.7      2020-11-20 [1] CRAN (R 4.0.2)
 pkgconfig     2.0.3      2019-09-22 [1] CRAN (R 4.0.0)
 plyr          1.8.6      2020-03-03 [1] CRAN (R 4.0.0)
 pROC          1.17.0.1   2021-01-13 [1] CRAN (R 4.0.2)
 prodlim       2019.11.13 2019-11-17 [1] CRAN (R 4.0.0)
 purrr       * 0.3.4      2020-04-17 [1] CRAN (R 4.0.0)
 R6            2.5.0      2020-10-28 [1] CRAN (R 4.0.2)
 Rcpp          1.0.6      2021-01-15 [1] CRAN (R 4.0.2)
 readr       * 1.4.0      2020-10-05 [1] CRAN (R 4.0.2)
 readxl        1.3.1      2019-03-13 [1] CRAN (R 4.0.0)
 recipes     * 0.1.15     2020-11-11 [1] CRAN (R 4.0.2)
 reprex        1.0.0      2021-01-27 [1] CRAN (R 4.0.3)
 rlang         0.4.10     2020-12-30 [1] CRAN (R 4.0.2)
 rmarkdown     2.6        2020-12-14 [1] CRAN (R 4.0.2)
 rmdformats  * 1.0.1      2021-01-13 [1] CRAN (R 4.0.2)
 rpart         4.1-15     2019-04-12 [1] CRAN (R 4.0.3)
 rprojroot     2.0.2      2020-11-15 [1] CRAN (R 4.0.2)
 rsample     * 0.0.8      2020-09-23 [1] CRAN (R 4.0.2)
 rstudioapi    0.13       2020-11-12 [1] CRAN (R 4.0.2)
 rvest         0.3.6      2020-07-25 [1] CRAN (R 4.0.2)
 scales      * 1.1.1      2020-05-11 [1] CRAN (R 4.0.0)
 sessioninfo   1.1.1      2018-11-05 [1] CRAN (R 4.0.0)
 stringi       1.5.3      2020-09-09 [1] CRAN (R 4.0.2)
 stringr     * 1.4.0      2019-02-10 [1] CRAN (R 4.0.0)
 survival      3.2-7      2020-09-28 [1] CRAN (R 4.0.3)
 tibble      * 3.0.5      2021-01-15 [1] CRAN (R 4.0.2)
 tidymodels  * 0.1.2      2020-11-22 [1] CRAN (R 4.0.2)
 tidyr       * 1.1.2      2020-08-27 [1] CRAN (R 4.0.2)
 tidyselect    1.1.0      2020-05-11 [1] CRAN (R 4.0.0)
 tidyverse   * 1.3.0      2019-11-21 [1] CRAN (R 4.0.0)
 timeDate      3043.102   2018-02-21 [1] CRAN (R 4.0.0)
 tune        * 0.1.2      2020-11-17 [1] CRAN (R 4.0.2)
 vctrs         0.3.6      2020-12-17 [1] CRAN (R 4.0.2)
 withr         2.4.1      2021-01-26 [1] CRAN (R 4.0.3)
 workflows   * 0.2.1      2020-10-08 [1] CRAN (R 4.0.2)
 xfun          0.20       2021-01-06 [1] CRAN (R 4.0.2)
 xml2          1.3.2      2020-04-23 [1] CRAN (R 4.0.0)
 yaml          2.2.1      2020-02-01 [1] CRAN (R 4.0.0)
 yardstick   * 0.0.7      2020-07-13 [1] CRAN (R 4.0.2)

[1] /Library/Frameworks/R.framework/Versions/4.0/Resources/library